function hpercentile = Hazenpercentile(data,perc)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   COMPUTE THE HAZEN PERCENTILE OF A RECREATIONAL WATER QUALITY DATASET                           
%   ASSUMING THAT ALL THE DATA ARE 'RWQ STANDARD' DATA, THAT IS: SAMPLE 
%   DATA TO BE USED IN THE HAZEN PERCENTILE CALCULATION. THESE DATA ARE                 
%   SUPPOSED TO BE TAGGED IN THE DATASET. TO SPECIFICALLY COMPUTE A 5-YR                   
%   95TH HAZEN PERCENTILE VALUE, USE 5YEAR OF  WEEKLY DATA. 
% THE RESULTING VALUE CORRSPONDS TO THE ANNUAL GRADE OF THE SITE   
%
% INPUT     : - data: data array                                                     
%             - perc: percentile value 
%
% OUTPUT    : - hpercentile : Hazen perenctile value
%  
% EXAMPLE: load TestData.mat  %(load 5 yr of weekly Enterococci data)
%          grade = Hazenpercentile(testdata.weekly5yr(:),95) % calculate
%                                    %  the 5yr 95th hazen percentile value
%                                                    
%          >grade = 91.8
%
% Author: Benedicte Pasquer, IMOS/eMII (http://imos.org.au/)
% email: benedicte.pasquer@utas.edu.au
% May 2013                  
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

data(isnan(data)|data ==99999)=[]; % remove missing values
[sorteddata,ind] = sort(data);
n = length(data);
p = perc/100;

if n >= 1/(2*(1-p)) &&  n >= 1/(2*p)
    rhazen = 0.5 +p*n;
    rI = floor(rhazen);
    rf = rhazen - rI;
    hpercentile = sorteddata(rI) + rf*(sorteddata(rI+1) - sorteddata(rI));
end

